Checking if the instalation went alright

Installation

  • To use pyLLE, you need a working python 3 installation. The easiest way (and best way in my opinion) is to install the anaconda distribution. When installing please check the option to add python to the path, as it would allows you to directly type the python command into a temminal/cmd window.
  • When Python and Julia 1.1.1 is installed, in a terminal (Linux/Mac Os) or a cmd window (Windows) type:
pip install pyLLE
  • This will automatically install the different python script and the needed Julia package. However, it is strongly advised to check if the installation went fine.

Cheking the Python Installation

  • Open a python command prompt. You can simply type python in a terminal (Linux/Mac Os) or a cmd window (Windows), or you can launch spyder, jupyter or a ipython instance.
  • Inside python, try to import the different package that pyLLE requires to check if they were installed correctly:
In [6]:
import scipy
import plotly
import numpy
import matplotlib
import h5py
import prettytable
import matplotlib
import ipdb

If any errors show up due to a missing package please:

  • Open a terminal (Linux/MacOs) or cmd window (Windows). The following commands should not be type in a python instance but directly in the shell.
  • Install the missing package by typing:
pip install <missing package name>

Checking the Julia Installation

  • open a Julia shell
  • type the following commands:
using Pkg

Pkg.add("HDF5")
Pkg.update("HDF5")

Pkg.add("FFTW")
Pkg.update("FFTW")

Pkg.add("LinearAlgebra")
Pkg.update("LinearAlgebra")
  • You should not have any errors at that level usually.
  • If you have an issue with the installation of a Julia package, directly enquire to the faulty pacakge through their github page as they might have solution to install it (the odd that is happens is really small though)

pyLLE example - Ring Resonator comb simulation

Please note that if you are using this notebook, you need to make sure that you "trust" it to ensure the display of the figure. To do so, in the jupyter notebook, go to File → Trust this notebook. It should reload the notebook and the plotly figure will be able to be displayed now

Import and Setup

In a python shell, spyder, a script, or in jupyter notebook, start by importing the package.

In [7]:
import pyLLE

We now define the resonator parameters. Here we will use a file TestDispersion.csv which is made of two columns: the first one is the azimuthal mode order (integer), the second is their corresponding frequency of resonance in Hz.

In [8]:
res = {'R': 23e-6, # ring radius in meter
       'Qi': 1e6,  # Intrinsic Q factor
       'Qc': 1e6,  # Coupled Q factor
       'γ': 1.55,  # Non-linear coefficient at the pump frequency
       'dispfile': 'TestDispersion.csv', # frequency and corresponding azymuthal mode simulated previously
      } 

We now define the simulation parameters. Here we precise a linear detuning ramp of the pump from δω_init to δω_end relative to the pump mode angular frequency, mode closest to the defined pump frequency f_pmp. The simulation length Tscan is in unit of round trip, as it is more convenient in the Lugiato-Lefever formalism. It is important to notice that two parameters for the mode bandwidth have to be defined, μ_fit which determined the fit window of the raw data found in dispfile, and μ_sim which is the number of mode simulated in the LLE, hence could be larger than the fit mode through extrapolation

In [17]:
import numpy as np
sim = {'Pin': 150e-3, # Input power in Q
       'Tscan': 1e6,  # Length of the simulation in unit of round trip 
       'f_pmp': 191e12, # Pump Frequency 
       'δω_init': 2e9*2*np.pi, # Initial detuning of the pump in rad/s 
       'δω_end': -8e9*2*np.pi,  # End detunin of the pump in rad/s
       'μ_sim': [-74,170],  # azimuthal mode to simulate on the left and right side of the pump
       'μ_fit': [-71, 180], # azimuthal mode to fit the dispersion on the left and right side of the pump
        }

In both the resonator and simulation dictionaries, the parameters can be called through their greek letters or through their equivalent latin names (e.g. μ → mu or δω → deltaomega). A translator dictionnary is implemented to translates every greek entries (see. self._greek)


Let's initialize the class

In [18]:
solver = pyLLE.LLEsolver(sim=sim,
                       res=res)

Dispersion Analyse

To plot and retrieve all the data of the dispersion, the method self.Analyze has to be called, resulting in the plot of the integrated dispersion Dint

In [19]:
dispfig = solver.Analyze(plot=True,plottype='all')
-- Dispersion Analysis --
	Pump index: 72
	Center Pump: 190.523 THz
	FSR: 1006.36 GHz
	D2: 24.58 MHz

One can clearly see that because we simulate a larger window with the LLE (orange curve - LLE simulation) than the raw data (blue dots), we extrapolate outside of this region. One has to be carefull about this feature as ripple in the integrated dispersion can happen causing zero-crossings which are artefacts


A new attribute disp has been created which consists of a dictionary of the different value of the retrieve dispersion

In [20]:
solver.disp.keys()
Out[20]:
dict_keys(['freq', 'ng', 'neff', 'D', 'Dint', 'dphi'])

Temporal Solver

One can solver the full temporal Lugiato Lefever equation :

$t_R \frac{\partial E(t, \tau)}{\partial t} = - \left(\frac{\alpha'}{2} - i\delta_0 \right)E + i \cdot \mathrm{FT}^{-1}\left[ -t_R D_{int}(\omega) \cdot \mathrm{FT}\left[E(t, \tau)\right]\right] + \gamma|E|^2 E + \sqrt{\theta}E_{in}$

The solver implemented in this package is based on a Julia core called through python. Hence to interface both language in an easy way, a hdf5 file is created in a temporary location with all the necessary data to solve the above LLE.

Hence, to solver the LLE here, we first need to setup the file:

In [21]:
solver.Setup()
-- Solving standard LLE --
	Simulation Parameters
		R = 23.00 µm
		Qi = 1.00 M
		Qc = 1.00 M
		γ = 1.55 
	Simulation Parameters
		Pin = 150.00 mW
		Tscan = 1.00 x1e6 Round Trip
		f_pmp = 191.00 THz
		δω_init = 2.00 x2π GHz
		δω_end = -8.00 x2π GHz
		μ_sim = [-74.00,170.00] 
		μ_fit = [-71.00,180.00] 


And now we can solve the equation

In [22]:
solver.SolveTemporal()
----------------------------------------------------------------------
2019-07-19 14:34:03
Launching Julia....
Temp file can be found in: /var/folders/qh/mzl16lsd5dj4254_sx6ghdsc0000gn/T/tmpl6ugg58clog.log
Launching Julia: Done
Computing LLE [**************************************************] 100%

Simulation Time 00h:08min:24.1s
----------------------------------------------------------------------

It is important to note that a KeyboardInterupt will result in killing the Julia process


When done, we need to retrieve the data

In [23]:
solver.RetrieveData()

Post Process of Temporal Solver

Here is a list of different method that ease the plotting of the quantity of interest, such as comb spectra and envelop of the electric field regarding the fast time. Although this is quite convenient for most of applications, one could use the raw complex electric field within the resonator. The different quantities saved by the solver are can be retrieved in the dictionary solver.sol

In [24]:
solver.sol.keys()
Out[24]:
dict_keys(['u_probe', 'Em_probe', 'Ewg', 'ω', 'comb_power', 'detuning', 'freq', 'theta'])

Everything in the solver.sol dictionarry is in SI units, except for the comb power, which is in dBm.

  • solver.sol['u_probe'], is the complex electric field relative to the fast time within the cavity, hence the complex number, probed by the solver. It is a matrix 1000 x [number of µ solved]
  • solver.sol['Em_probe'] is the FFT of solver.sol['u_probe'], hence the complex spectra of the electric field. The comb spectra, within the cavity is given bu $10\times\mathrm{log10}(|\mathtt{Em\_probe}|^2)$ From Em_probe, it is easy to find the output comb power, given by Em_probe*sqrt(η), where η (in the Julia solver), is the coupling rate to the waveguide, given by Qc
  • solver.sol['Ewg'] corresponding to $\mathtt{solver.sol['Em_probe']}\times\sqrt{\eta}$, hence the complex spectra envelop of the electric field at the output of the coupled waveguide.
  • solver.sol['detuning'] is the detuning ramp used
  • solver.sol['ω'] is the angular frequency
  • solver.sol['ω'] is the optical frequency
  • solver.sol['theta'] is the position within the cavity, in radian (between -π and +π).

It is import to note that all these quantity are subsampled in regards of the split-Fourier total step taken. Indeed, the number of step to solver the temporal LLE would be to large to efficienlty save or display the results. Hence, we subsample the steps taken by the solver at a fix number of 1000 (hence why the different matrices have a size of 1000x [number of µ solved])

Displaying Results

We can plot the first figure that gives an overall idea of the behavior of the resonator

In [25]:
fig_sol = solver.PlotCombPower()

To have a better idea of what happen for a given detuning we can plot the spectra and the temporal profile of the electric field in the resonator

In [27]:
ind = 704
fig_temp = solver.PlotCombSpectra(ind)

The temporal profile can also be retrieve, interesting in a case of we are on a soliton state for the given detuning

In [28]:
_ = solver.PlotSolitonTime(ind)

Saving figures

Helper to save the figures, in any format supported by matplotlib have been created. The helper (_pylle._llesolver_Latexify_) has been designed to make "publication-ready-like" figure, while one is using jupyter notebook/jupyter lab (hence reliying on plotly for figure display) or through a python consol (hence thorugh matplotlib display). Depending on which type of plot has already been display (comb summary, comb spectra, time profile), this figures will be saved in the specified format

In [ ]:
solver.SavePlots2File(basename = 'images/', format = 'pdf')
solver.SavePlots2File(basename = 'images/', format = 'png')

Pickling the solver

The whole package has been design with simplicity in mind, which mean sinplicity to compute but also simplicity to retrieve results. Hence pyLLE is coded in a way that, using the pickle package, one can easily save the state of the solver and retrieve it later, with all the feature presented above still available. It is important to use the self.SaveResults method to pickle the solver as it takes care of dumping the different thread that cannot be pickle in the class.

In [29]:
solver.SaveResults('PickleSolver', path = './')
./PickleSolver.pkl

One can reload the state of the simualtion:

In [30]:
import pickle as pkl
oldSolver = pkl.load(open('PickleSolver.pkl', 'br'))

The different methods can be called with the loaded pickled solver

In [31]:
figly = oldSolver.PlotCombSpectra(740)

Steady state Solver

One can solver the steady state Lugiato Lefever equation :

$- \left(\frac{\alpha'}{2} - i\delta_0 \right)E + i \cdot \mathrm{FT}^{-1}\left[ -t_R D_{int}(\omega) \cdot \mathrm{FT}\left[E(t, \tau)\right]\right] + \gamma|E|^2 E + \sqrt{\theta}E_{in} = 0$

Using a Newton method to find the root of the equation will result in solvin this particular state of the LLE.

The method self.SolveSteadySteate is design for this.

Although, it gives fast results, the accuracy of such solver remains questionable compare to a full temporal resolution of the equation


First we need to change the simulation parameters to introduce a fixe detuning δω

In [32]:
solver.sim['δω'] = -3.6e9*2*np.pi # more or less what it is as the end of the soliton step 
In [34]:
steady_fig = solver.SolveSteadyState()
----------------------------------------------------------------------
2019-07-19 14:49:18
Not symmetric mode calculation -> switching to it with µ_sim = [-170, 170]

-- Dispersion Analysis --
	Pump index: 72
	Center Pump: 190.523 THz
	FSR: 1006.36 GHz
	D2: 24.58 MHz